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Abstract 

We consider the phenomenon of spectral pollution arising in calcu- 
lation of spectra of self-adjoint operators by projection methods. We 
suggest a strategy of dealing with spectral pollution by using the so- 
called second order relative spectra. The effectiveness of the method 
is illustrated by a detailed analysis of two model examples. 
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1 Introduction 



Let A be a self-adjoint operator on a Hilbert space H. In order to de- 
termine the spectrum Spec (A) of A, one can use the following projection 
method. Let (Ck)ke'N be a sequence of closed linear subspaces of the domain 
Dom(A) of A and suppose that the corresponding orthogonal projections 
Pk '■ Ti. — > C-k converge strongly to the identity operator /. Let A (A) be the 
set of all such sequences of subspaces. One may ask whether or not 

lim Spec(v4, = Spec(A), (1) 

where Spec(A, Cu) is the spectrum of P^A : Ck —* C-k and "lim" is defined 
in an appropriate sense. Unfortunately the answer is, in general, negative 
and lim Spec(v4,£fc) may contain points which do not belong to Spec(A). 

This phenomenon is well known in numerical analysis and is called spectral 
pollution. 

Since we allow A to be unbounded, it is convenient for us to extend 
Spec(A) by adding to it — oo, +00, or both, if A is unbounded below, above 
or from both sides, respectively. We denote the resulting set by Spec(A). In 
other words, Spec (A) is the closure of Spec (A) in the two point compactifi- 
cation 1 := RU {±00} of BL Let Spec ess (v4) be the essential spectrum of A, 
i.e. 

Spec css (A) := Spec(A) \ {isolated eigenvalues of finite multiplicity}, 

and let Spec ess (A) be the extended essential spectrum of A, i.e. 

Spec ess (A) := Spec(v4) \ {isolated eigenvalues of finite multiplicity}. 

It turns out that spectral pollution may occur at any point lying in 
a gap of the extended essential spectrum, i.e. at any point of the set 
conv ^Spec ess (A)j \ Spec css (A), where "conv" means the convex hull (see 

the proof of Theorem [O] below, where this fact is demonstrated by choosing 
an appropriate sequence of finite-dimensional subspaces) . 



One can argue that the sequence (Ck)keN m the proof of Theorem |2J] had 
been deliberately chosen in a particularly nasty way and that things like this 
do not usually happen in "real life". This argument is only partially true 
at best, and examples in Sections 3 and 4 show how the standard approach 
leads to spectral pollution in some very straightforward problems. 

One can try to deal with spectral pollution by choosing C k in a "reason- 
able way". However the known recipes (see, e.g., [ |BBO[ ], | |BDO| |, [ |RSHSP V 
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and references therein) are by no means universal or rigorously justified for a 
wide range of problems. An absolutely "safe" choice of (Ck)km would be the 
one based on the spectral decomposition of A. Unfortunately knowing the 
spectral decomposition is much more than knowing just the spectrum which 
is to be found. 

We propose a complementary (not an alternative!) way of dealing with 
spectral pollution, which is based on the notion of the second order spectrum 
Spec 2 (A,£) of A relative to a subspace C C Dom(A) (see Definition f2~3| ) 



introduced in |Davj| . Second order relative spectra usually contain complex 



numbers, but it turns out (see Theorem |2.5|) that if z G Spec 2 (A, C) then 



Spec(A) n [Rez- | lmz\,Rez + | Imz|] ^ 0. (2) 

In particular, this means that if z belongs to a second order relative spectrum 
of A, and its imaginary part is small, than there is a point of the spectrum 
of A close to Re z. 



Q 



Re z-llm z| 



Re z+llm z| 



Figure 1: Spec(^) n [Rez - | Imz|,Rez + | Imz|] / 0. 

This fact suggests the following strategy, which uses @ as a posteriori 
estimate. One should choose a reasonable sequence (Ck)ken and apply the 
standard projection method, i.e. find Spec(v4, £&). Additionally, one should 
find Spec 2 (A Ck) (for possibly different Ck C Dom(A) ). Suppose Spec(A, Ck) 
contain points which appear to converge to some A G M as k — > oo. // 
Spec 2 (A, Ck) also contain points lying close to X, then X is indeed close to the 
spectrum of A. If however Spec 2 (A, Ck) do not contain points approximating 
X, one cannot conclude that X is far from the spectrum. A definite claim 
which can be made in this case is that X requires additional attention^ as 

1 One can, for exa mpl e, use the method suggested in the companion paper by E.B. 
Davies and M. Plum fDFj. 
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spectral pollution may be around the corner. It has to be said though, that 
our (very limited!) experience seems to suggest that Spec 2 (A, £ k ) do capture 
the spectrum quite well. 

We do not claim that the above method is going to be universally suc- 
cessful. Our main ambition is merely to bring it to the attention of those 
concerned with finding spectra numerically. We have applied it to two prob- 
lems which were chosen mainly for their simplicity. The results (which turned 
out to be better than we expected) are presented in Sections 3, 4, 6, and 7. 

Section 2 contains the main theoretical results on standard (first order) 
and second order relative spectra. Similar results on second (and higher) 
order relative spectra have been obtained in ||Sha|| under slightly stronger 
restrictions. We start with Theorem |2.1| , which has been mentioned above 
and which says that spectral pollution may occur at any point lying in a gap 
of the extended essential spectrum. It implies that one should be particu- 
larly careful when applying a projection method to find the spectrum of an 
operator A in the following cases: 

a) A is not semi-bounded and Spec (A) ^ R; 

b) A is semi-bounded but unbounded, and has a nonempty essential spec- 
trum different from a half-line; 

c) A is bounded and has at least one gap in the essential spectrum. 
Examples of types c) and b) are treated in Sections 3 and 4 respectively. We 
return to these examples in Sections 6 and 7 after giving the proofs of the 
main theoretical results in Section 5. 

Section 2 also discusses an additional difficulty which arises for unbounded 
operators. It is well known that if A is bounded, then the left-hand side of 

d) ]) contains the right-hand side. If, however, A is unbounded above, say, 
then the left-hand side of (JlJ) may well be equal to {+00} (see Theorem |2.2|) . 

In order to make the paper more readable, we have divided the material 
into two parts. The first one consists of Sections 2-4 and covers the main 
theoretical and numerical results. The same material is treated at a greater 
depth in the second part of the paper, which consists of Sections 5-7. The 
summary is given in Section 8. 

Finally, the accompanying webpage 

www .ma . hw . ac . uk/ ~levitin/SpectralPollution 
contains the listings of the Matlab programmes written by the authors and 
used for numerical calculations in Examples I and II. 
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2 Main results 

Let A be a self-adjoint operator on a Hilbert space 7i. Throughout the 
paper E(-) denote the spectral projections associated with A. The proofs of 
all results of this Section are given in Section 5. 

The following result states that spectral pollution may occur at any point 
in a gap of the extended essential spectrum of A. 

Theorem 2.1. For any A G conv ^Spec ess (A)j \ Spec css (A) there exists an 
increasing sequence (Ck)keN £ M^.) suc ^ ^at 

AeSpec(AA), VfceN. 

It is well known that for bounded self-adjoint operators the left-hand side 
of (0) contains the right-hand side: 



lim Spec(A, Ck) ^ Spec(A) (3) 



see, e.g., |[Arv|| or [ |Sha| ] ) . So, spectral pollution means that the left-hand side 



of (0) is strictly larger than the right-hand side. For unbounded self-adjoint 
operators the situation is more complicated: the left-hand side of ([!]) does not 
necessarily contain the right-hand side. The reason for this is intuitively very 
clear. Let be finite dimensional and suppose for definiteness that A is not 
bounded above. Then perturbing a basis of by arbitrarily small vectors 
lying in Dom(v4) fl E([N, +oo))7i with sufficiently large N, one can push the 
spectrum of the corresponding finite dimensional operator arbitrarily far to 
the right. 

Theorem 2.2. Let A be an unbounded above self-adjoint operator on a sep- 
arable Hilbert space TC. Then for an arbitrary sequence (£fc)fceN G A(A) 
of finite dimensional subspaces and arbitrary sequences (ek)ken \ and 
(Rk)keN +oo, there exists a sequence (C k )k£N £ -M^4) suc h that 

\\P k - P' k \\ < s k and Spec(A,^) C (R k ,+oo), VA; 6 N, 



where P' k : TC — * C' k are the corresponding orthogonal projections. A similar 
statement holds for operators unbounded below. 

The requirement C Dom(A) is too restrictive for many applications, in 
particular for approximating spectra of differential operators by finite element 
methods. When using finite element methods, one usually states the problem 
in a weak (variational) form and requires that belong to the domain of 
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the corresponding sesquilinear form. The restriction Ck C Dom(A) would 
mean extra smoothness of the finite elements. 

Suppose A is a bounded below self-adjoint operator on a separable Hilbert 
space 7i and let Qa '■ 7i x TC — ► C be the corresponding closed sesquilinear 



form (see, e.g., ||Kat| , §1, Ch. VI]). Let (Ck)keN be a sequence of finite di- 
mensional linear subspaces of the domain Dom(Q J 4) of Qa such that the 
corresponding orthogonal projections Pk : H — > Ck converge strongly to /. 
Let A(Qa) be the set of all such sequences. A number A G C is said to belong 
to Spec(A,£ fc ) if there exists u G £& \ {0} such that 

Qa(u, v) - X(u, v) = 0, \/v G £ fe . 

If Ck C Dom(y4), then Qa{u,v) = (Au,v) and the above definition means 
that 

• there exists u G £& \ {0} such that (Au, v) — X(u, v) — 0, \fv G £fc, 

• i.e. there exists u E Ck \ {0} such that P^Azz = Am 

• i.e. A belongs to the spectrum of PkA : Ck — > Ck- 

Hence one gets the old definition of Spec(A, Ck) for finite dimensional Ck C 
Dom(A) (cf. Remarks ^.4| and |2.6|(b) below). 



It is not difficult to prove that Theorem |27| remains true if A (A) is re- 
placed by A(Qa) (see Section 5). 

The above results show the difficulties one may face when trying to ap- 
proximate the spectrum of a self-adjoint operator with the help of a projec- 
tion method. A possible way of addressing spectral pollution is based on the 
notion of the second order relative spectrum. 

Definition 2.3. Let A be a self-adjoint operator on a Hilbert space H and 
let C be a finite dimensional subspace of Dom(A) . A number z G C is said to 
belong to the second order spectrum Spec 2 (^4, C) of A relative to C if there 
exists u G C \ {0} such that 

{{A — zl)u, {A — zl)v) — 0, Wg£, (4) 

where (*, •) is the scalar product in 7i. 

Remark 2.4. Suppose C C Dom(A 2 ) and let P be the orthogonal projection 
onto C. Then z G Spec 2 (A,£) 

• iff there exists u G £ \ {0} such that ((A — zl) 2 u, v) = 0, Vt> G £, 

• i.e. iff there exists u G £ \ {0} such that u G Kei(P(A — zl) 2 ), 

• i.e. iff the operator P{A — zl) 2 : C — > C is not invertible. 

The last is the definition of the second order relative spectrum introduced in 
Davfl and used in IS ha . 
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Theorem 2.5. If z E Spec 2 (A, £) then 



Spec(A) n [Re z- | Im;z|,Rez+ | hnz\] ^ 0. 



Remark 2.6. (a) Theorem |2.5| has been proved in [ 5h"a in the case C C 
Dom(A 2 ). 

(b) Everything said above about second order relative spectra remains un- 
changed if £ is infinite dimensional, with only one exception: the implication 



P(A — zl) 2 : C — > C is not invertible 
(see Remark |2.4|) is no longer true. 



Ker(P(A - zlf) n C + {0} 



3 Example I: operator of multiplication by a 
discontinuous function 

Consider the operator of multiplication by a function a: 
A = al : L 2 ([-7r, it]) -> L 2 ([-7r, vr]) , 

where 

( 3 1 / - 

, . I h - cos v 5 x , for — 7r < a; < , 

a(x) = < 2 2 

2 + cos V^x, for 0< 

It is clear that 

Spec(A) = Spec css (A) = a ([-tt, it}) = [-2, -1] U [1, 3] . 

The function a and, in particular, the irrational factors a/5 and y/2 have been 
chosen to get a "generic" situation and avoid any possible hidden symmetries. 

Using the standard orthonormal basis e^{x) := (27r)~ 1//2 exp (ikx), k e Z, 
one can represent A as a (doubly) infinite Toeplitz matrix and take the 
orthogonal projections onto Cjy := span ({efc}£ = _„), where iV = dim/^ = 
2n + 1. Figure ^| shows the approximation of the spectrum obtained by 
the standard projection method for a typical value of N. One can see that 
the essential spectrum [—2, —1] U [1,3] is approximated very well, but that, 
additionally, there are some computed eigenvalues lying in the gap of the 
essential spectrum. 

The spurious eigenvalues vary very slowly with N, see Table [I], and their 
number #(Spec(A, C N ) n (—1,1)) also increases slowly as iV increases. In 
fact, they seem to converge monotonically to some points in (—1,1), and 
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Figure 2: Spec(AAooi)- 



N 


101 401 701 1001 


Spec(A,£ N ) n (-1,1) 


-0.9997 -0.9985 
-0.9834 -0.9673 -0.9550 
-0.9264 -0.8250 -0.7823 -0.7553 
-0.4866 -0.3478 -0.3002 -0.2721 
0.4362 0.4597 0.4673 0.4717 



Table 1: Eigenvalues lying in the gap of the essential spectrum of A for some 
values of N. 



judging by the numerical evidence alone, one might have thought that the 
limits belonged to the spectrum, if the example was less trivial and the 
spectrum was not known. 

The presence of spurious points is hardly surprising. Indeed, if instead of 
the operator of multiplication by a, one considers the corresponding Toeplitz 
operator T a acting on the Hardy space H 2 = span ({ek}k>o), then the projec- 
tion method using A4n '■= span ({e^}^ 1 ) leads to the same finite-dimensional 
Toeplitz matrices. In other words, the standard projection method does not 
know which spectrum is being calculated: that of the operator of multipli- 
cation or of the Toeplitz operator. The latter equals [—2,3] (see [tiW| or 



|[BS|, 2.36]). Therefore not only will we have spectral pollution, the spurious 



eigenvalues will eventually fill the entire gap (—1,1) (cf. (j3|) or ||Bot|| ). 

Remark 3.1. It is interesting that spurious eigenvalues appear so slowly. 
Three spurious eigenvalues are present already for rather small values of N. 
The fourth one emerges from —1 for iV pa 161, and the fifth one emerges from 
— 1 for N m 631. Such a slow rate seems to have the following explanation. 
The part (—1,1) of the spectrum of the Toeplitz operator T a is there because 
of the discontinuities of a(x) at x = and x = —tt ~ tt. Functions responsible 
for this part of the spectrum appear to be very localized at and — 7r ~ 7r 
and hence have a very slow decay of the Fourier coefficients. Figure |3] shows 
the graph of the trigonometric polynomial which is an eigenfunction of PnA : 
C N -> C N , N = 1001, corresponding to A w 0.4717. 
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1 



o - 



-1 - 



-2 - 



-3 - 



-4 - 



-5 - 




-4-3-2-1 1 2 3 4 

Figure 3: A typical computed "eigenfunction" of A corresponding to a spurious 
eigenvalue. 

We now show how second order relative spectra locate the spectrum of 
A and, in particular, detect spectral pollution in the standard projection 
method. 

Figure f| shows the second order relative spectra Spec 2 (^4, £n) computed 
for iV = 101,201,501. 

For each N, the second order relative spectrum consists of a few points 
lying approximately on a circle centred on the real axis and two contours "sur- 
rounding" two intervals of the (essential) spectrum. The number of points 
of Spec 2 (A C N ) lying on the circle does not seem to increase with N. 

Figure |^ gives a detailed view of Figure |] near the real axis (note the 
different scales along the coordinate axes). The contours made of the points 
of the second order relative spectra, which "surround" two intervals of the 
(essential) spectrum, get closer to the real line as N increases. Most im- 
portantly, the absence of points of the second order relative spectra close to 
the eigenvalues found by the standard projection method inside the interval 
(—1,1) suggests that these eigenvalues are spurious. 
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Figure 4: Spec 2 (A, £ N ). 



In order to illustrate a situation when an operator has both non-empty 
essential and discrete spectra, let us consider a two-dimensional perturbation 
of A = aJ, 



Bu(x) := a(x)u(x) / u(t) exp (it)dt 

27T 



+ expj>x) f u ^ exp (_ 2 ^) t ft j xe[— 7T,7r] 
27T y_ 7r 



(5) 



This time, the essential spectrum is again very well approximated by the 
usual projection method, and the perturbation gives rise to several eigenval- 
ues, see Figure |], listed here for N = 401 : one (Ai ~ —3.6057) to the left 
of the essential spectrum, six (A 2 « -0.9842, A 3 « -0.8947, A 4 « -0.8203, 
A 5 « -0.3496, A 6 « 0.4603, and A 7 « 0.6082) in the gap (-1,1), and one 
Ag ~ 3.0432) to the right of the essential spectrum. 

To establish which of the eight eigenvalues are spurious, we compute the 
second order relative spectra of B. Their structure is similar to that of 
Spec 2 (A,£7v) apart from several additional points located very close to the 
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Figure 5: Detailed view of Figure near the real axis. 



N= 401 



Figure 6: Spec(S, £401). 



real axis, see Figures |7| and §]. In particular, for N = 401 these points are 
approximately -3.6056 ± 0.0505i, -0.8929 ± 0.0266i, 0.6085 ± 0.0442i, and 
3.0413±0.0178i. Thus, the eigenvalues Ai, A3, A 7 and A 8 are correct, whereas 
the other eigenvalues A2, A4, A5, and A6 are almost certainly spurious. 
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Figure 7: Spec 2 (B,C 
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= 401 


• Spec(B) for N = 


401 




-0.15 - 



-0.2 1 1 1 1 1 1 1 1 
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Figure 8: Detailed view of Figure f?] near the real axis and, for comparison, 
Spec(l?,£Ar). 
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4 Example II: Stokes type system of ODEs 



Our next example is an operator H in L 2 ([— tt, 7r]) generated by a simple 
2x2 Stokes type system of ODEs subject to the periodic boundary conditions, 



/ d 2 



H 



dx 2 



d \ 

dx 



, ± 2 

V dx 



J 



Finding the spectrum of H explicitly is easy. The essential spectrum is found 
by using a general result on Agmon-Douglis-Nirenberg (that is, mixed order) 
elliptic systems ||GG|| , which gives 



Spec ess (#) = AGC 3(Gf\{0}: det 



e 



-i£ 2 - A 

To find the discrete spectrum, consider the system of equations 



{I}- 



-u" - v 



Xu 

u' + 2v = \v 



(6) 



subject to the periodic boundary conditions. Expressing v from the second 
equation and substituting the result into the first one gives 



—u 



A(A - 2) 
A- 1 



u 



(7) 



Comparing this equation with the standard periodic spectral problem for 
—d 2 /dx 2 we conclude that 

A G Spec disc (#) iff X ^ ~ ^ G Spec ( = {k 2 : k E Z} . (8) 



A- 1 



This implies that H has two series of eigenvalues A^ given by 



A 



± 



k 2 + 2 ± ^(k 2 + 2) 2 - 4fc 2 



k = 0,1,2, 



Each eigenvalue A fc or A^ with k > is of multiplicity 2. The eigenvalues 
of the "minus" series are located in the interval [0, 1) and converge to 1 
as k — > oo. The eigenvalues of the "plus" series are located in the interval 
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[2, +oo) and tend to +00 as k — > 00. Note that the interval (1, 2) is spectrum- 



free. The eigenvectors (u,v) of H corresponding to the eigenvalue A fc , k > 



are linear combinations of 
1 

<fc/(Af-2)r \^k/(X 



1 ^d ( . * Je- fe . (9) 



A' 

+ 



The eigenvectors corresponding to A =0 and \q = 2 are proportional to 
* and M respectively. 



°/ V 1 , 

The projection method we use in the numerical analysis of the spectrum 
of H is the standard finite element method (FEM). It turns out (see Section 
7 below) that the results of calculations depend dramatically upon the choice 
of the corresponding spaces of finite elements. 

We start with the simplest possible choice — the space of one dimen- 
sional periodic Lagrange elements of order p, that is, of continuous functions 
which are piecewise polynomials of degree p. Here N is the number of mesh 
node points on the interval [— ir, ir], so the total dimension of is Np. 

Remark 4.1. In most of the numerical experiments, we use just the linear 
Lagrange elements (i.e., of order 1). We shall therefore write := Ljy. Ad- 
ditional tests have shown that varying the order affects neither the qualitative 
picture, nor the rate of convergence. 

It is important to note that the spaces If N are embedded into the Sobolev 
space VT 2 1 ([— 7r, 7r]) but not into the Sobolev spaces H /7 2 2 ([— 7r, 7r]). Thus, the 
elements of Ltv X Ltv do not necessarily belong to Dom(if). However, they 
belong to the domain of the corresponding quadratic form, which is sufficient 
for the FEM. 

The Lagrange elements are provided by all standard FEM packages, in- 
cluding, e.g., FEMLAB (see flFEMU ). 

The calculation of Spec(if, x Ljy) produces a very striking result: it 
gives very accurately several initial eigenvalues from the A 1 * 1 series, filling 
at the same time the entire interval (1, 2) with spurious eigenvalues (see 
Figure || which, for brevity, shows only the eigenvalues lying in the interval 
[0,6]). In the particular case N = 100, there are 15 eigenvalues in [0, 1], 99 
eigenvalues in [2, 00), and the remaining 86 eigenvalues are in (1, 2). 

The subspace £ = ~L N x L N used above is unsuitable for the calculation 
of the second order relative spectra Spec 2 (H, C) since for elements of this 
space the left-hand side of (f|) does not make sense. Therefore, as usual 
when dealing with fourth order ODES, let us introduce the standard space 
of periodic cubic Hermite elements M^, whose functions are continuously 
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Figure 9: FEM calculated Spec(i?, Ljy x Ljv) plotted against N. All eigenvalues 
in the interval (1,2) are spurious! 



differentiable and are piecewise polynomials of degree 3. They belong to the 
Sobolev space VFf ([— 7r, 7t]) so they are smoother than elements of Ljv- Here 
again iV is the number of mesh node points, so that dim Mat = 2N. 

Since the operator H is of mixed order, it is sufficient to choose the 
space Mat x Ljv for the projection space in (f|). The results of calculation of 
Spec 2 (H, Mjv x Ljv) are shown in Figures 10 and |ll[ The first figure (with 
different scales along the real and imaginary axes) demonstrates once more 
the advantage of this method in eliminating the spectral pollution. Note 
also the series of points of Spec 2 (iJ, Mjy x L N ) along the line Re A = 1 - 
this indicates the presence of the essential spectrum of H at A = 1. The 
second figure shows all computed points of the second order relative spectra 
for N = 25, 50, and 100. 

The points of the second order relative spectra give good approximations 
for the low eigenvalues of H. We postpone the discussion of the rates of 
convergence in various cases until Section 7. We note only that the number 
of well approximated eigenvalues increases with N, see Figure [TT| 
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Spec 2 (H,M N xL N )for/V=50 
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Figure 10: Second order relative spectra Spec 2 (-ff, Mat x Lat) for various values 
of N — low eigenvalues. 

Examples I and II are discussed further in Sections 6 and 7; we proceed 
now to the proofs of the main theoretical results stated in Section 2. 
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Figure 11: Second order relative spectra Spec 2 (-ff, Mjv x L^r) for various values 
of iV — all eigenvalues. 
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5 Proofs of the main results 



Proof of Theorem It is clear that there exists A_ G Spec ess (A) which lies 
strictly to the left of A. Consider first the case A_ G Spec css (A). If A__ is 
an isolated point of Spec(v4), then it is an eigenvalue of infinite multiplicity 
and we define Hj := span({ej}), where ej, j G N are mutually orthogonal 
eigenvectors of A corresponding to A_. If A_ is not an isolated point, there 
exist pair-wise disjoint compact sets W~ C M lying strictly to the left of A and 
such that Wj~ nSpec(A) /0,jG N. In this case we choose Hj = E(Wf)H, 
where E(-) are the spectral projections corresponding to A. Finally, if no 
point of Spec ess (^4) lies strictly to the left of A, then A_ = — oo and there 
exist eigenvalues Xj < A of A converging to — oo. Here we set Hj := E({Xj}). 
In all the above cases Hj C E((—oo, A)), j G N. 

Similarly one can define subspaces Hf C E((X, +00)), j G N. It is clear 
that Hf are all mutually orthogonal, Hf C Dom(A), AHf C Hf, and there 
exist r k > such that r& /* +00 and Hf C E([— r^, Vk])H, j < k. Let Ho be 
the closed linear subspace spanned by Hf, j G N and Hq be its orthogonal 
complement. 

Let vf G Hf, \\vf\\ = 1 and fif := (Avf,vf). Then it follows from the 
spectral theorem that n k < A < nf, k G N. Let Ck be the closed linear 
subspace spanned by Hf, j < k, E([— rk, t^Hq , and the vector 



„ • I Vt- X „- 1 / A - »k „+ 

Uk — \l — ~ v k + \ — v k . 

It is easily seen that the sequence {Lk)k&\ is increasing and belongs to A(A). 
Since vf,Avf G Hf C E([— rk, Tk])H, the orthogonality properties of Hf 
imply 

IKII = i> 

((A - XI)uk,u k ) = ^-^{Av k ,v- k ) + -^f-^(Avf,vf) - 

^k 



f4 










A- 


v k ■+ 






f4 - 



AIM 2 = -^i—^k + -r^fF^ A = ( J- 

{{A - XI)u k ,v) = 0, Vf G £([-r*,r fc ])?# U (J(WJ U Hf). 

j<k 

Hence ((A—XI)uk,u) = 0, Wu G i.e. Pfc(A— XI)uk = 0, i.e. PkAuk = Xuk, 
i.e. A is an eigenvalue of P^A : C k — > Vfc G N. □ 



Remark 5.1. If H is separable, the subspaces C k in Theorem |2J] can be 



made finite dimensional. Indeed, let e^ m , m = 1, . . . be an orthonormal 
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basis of Tif and e^, m — 1, . . . be an orthonormal basis of H.q. These bases 
exist because all subspaces of H are separable. Let Ck be the closed linear 
subspace spanned by ef m and E{[— r k , r fc ])e^ with j, m < k, and the vector 



\4 ~ A 



e k,i + 



Then it is clear that £& is finite dimensional, and the rest of the proof of 
Theorem |2.1| remains unchanged. 



Proof of Theorem \2.Q . It is sufficient to prove that for any finite dimensional 
£ C Dom(yl) and any e, R > there exists C C Dom(v4) such that 

||P - P'\\ < £ and Spec(A, £) C {R, +00), 

where P : 7i — > £ and P' : li. —> CJ are the corresponding orthogonal 
projections. 

Let m := dim£ and </>!,..., </? m be an orthonormal basis of C Since A 
is not bounded above, for any r > there exist mutually orthogonal vectors 
■01,..., ij) m G Dom(A) n E((t, +oo))TC which we normalize by ||^|| = 5, 
I = 1, . . . , m, where 5 G (0, 1). The numbers r and 5 will be chosen later. 

It is easy to see that the vectors <pi + fa, . . . , <p m + ip m are linearly inde- 
pendent. Indeed, for any c\, . . . , c m G C Parseval's identity implies 



^2ci((pi+fa) 



1=1 



> 



in 

^2 cm 

1=1 



E Q ^ 



1=1 



\ 1/2 

2 



(10) 



. 2=1 



Similarly 



^2ci(ipi+fa) 



1=1 



<(! + *) E 



1/2 



Q 



:n) 



j=i 



Let £' be the subspace spanned by <£>i + fa, . . . , yj m + ^ and consider 
the operators J : C — > £' and </':£'—►£ defined by 



J E c ^ = E c '(^ j + J/ E c ^ + ^) = E 



,«=i / ;=i 

It is clear that J'J = Jr, J J' = I a an d 



J=i 



(/ - j) Ij2 w = - E c ^ ( J - J ') E + ^) = E 



1=1 



1=1 
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It follows from Parseval's identity and (flQ|) , ([11]) that 

\\J\\<l + 6, ||J'||<— L- , \\{I-J)P\\<8, ||(/- J')P'\\ < * 



1-5' llv ' 11 " ' llv y 11 " 1-6 
Further, 

\\p - p'\\ < ||p - p'p\\ + \\p'p - p'\\ = ||p - p'p\\ + ||pp' - p'||, 

since (P'P — P')* = PP' — P' . Further, for any x G TC the closest point of £' 
isP'x. Hence ||P-P'P|| < ||P-JP||. Similarly \\P'-PP'\\ < \\P'-J'P'\\. 
Consequently 

hp - p'h < - j)p\\ + - j')p>\\ <t+Y^- 5 = 6 ~y^t- 

Fix b G (0, 1) such that 5(2 - 5)/{l - b) < e. Then ||P - P'|| < e. 

Since £ is a finite dimensional subspace of Dom(A), there exists M > 
such that ||Au|| < M||w||, Wu G C Take an arbitrary 



fit 

v 

i=i 



and set 

m m 

p:=Yl C M G A V> : = X] Q ^ G Dom(A) n P((r, +oo))W. 
i=i i=i 

Then 

(Au, v) = (Aip, ip)+2 Re(Acp, ip) + (Acp, <p) > 



r||^|| 2 - 2M\\(p\\ - M\\<p\\ 2 = (rb 2 - 2Mb - M) 



|2 



«=1 



Choose r such that (r<5 2 - 2Mb - M)/(l + b) 2 >R+l. Then 
/ „ x r<5 2 — 2Mb — M ll2 . ... .., 

^ — (TTbJ 2 — IHI ~ {R+ 1)IHI 

(see (|TTD). Hence Spec(A, C) C [R + 1, +oo). □ 

Suppose A is a bounded below self-adjoint operator. Then everything 
said above generalizes from A(A) to K(Qa) (see Section 2 for the definition) 
with almost no changes. One only needs to modify the concluding estimates 
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of the proof of Theorem |2.2| in the following way. Adding if necessary a 
sufficiently large positive number to the bounded below operator A, one can 
assume without loss of generality that A is positive definite. Then 

Q A (v, v) = Q A (ip, 1>) + 2 Re Q A (<p, ip) + Qa(<P, v) > 
Qa(M) -2Q a (V^) 1/2 Qa(<^) 1/2 > 

\q a (?i>, 40 - 2Qa(v, <p)>Iw\\ 2 - 2Myf = 

t6 2 -4M '"' 



E 



1=1 

Choose r such that (t5 2 - 4M)/2(1 +6) 2 >R+1. Then 

OiM > 2(1 + <y)2 IN > (* + l)|M| a . 

Hence Spec(A, £') C [i2 + 1, +oo). 

Now we proceed to the proof of Theorem [2.5| . For any c\,C2 G C we 
denote by D(ci,c 2 ) the open disk with diameter (ci,c 2 ). 

Lemma 5.2. Let A be a self-adjoint operator on a Hilbert space 7i and let L 
be a finite dimensional subspace o/Dom(A). Suppose a, b G M are such that 

(a, b) H Spec(A) = 0. 

Then 

D(o,6) n Spec 2 (A,£) = 0- 

Proof of Lemma \5. 4 Take an arbitrary 2; G D(a, 6). The interval (a, b) is 
seen from z at an angle ir/2 + e, where e G (0, vr/2]. It is clear that the set 

F = {(A - z) 2 : AG Spec(A)} 

lies in a sector with vertex at and angle it — 2e. The distance from to F 
is obviously positive. Therefore there exist 1? G (— tt, 7r] and c > such that 

Re (e**(A - z) 2 ) > c, VA G Spec(A). 
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Now applying the spectral theorem for A one obtains 

Re [e li) {{A - zl)u, {A - zl)u) ) = 

Re (e v& (||Au|| 2 - 2z(Au, u) + z 2 ||u|| 2 ) ) = 

Re ( e l ® [ (A 2 - 2Xz + z 2 )d(E(X)u, u)) = 

Re (e^(\ - zf) d{E{X)u, u) > c I d{E{X)u, u) = 

Spec(A) Jspec^) 

c||w|| 2 , Vw G C. 

Hence u G C \ {0} cannot satisfy (||), i.e. z G" Spec 2 (A, C). □ 

Proof of Theorem Suppose the intersection in the statement of the The- 
orem is empty. Then there exists e > such that 

Spec(A) fl [Re z — \ lmz\ — e, Rez + \ lmz\ + e] = 0. 

Hence D(Re z — | Im z\ — e, Re z + \ Im z\ + e) does not intersect Spec 2 (A C) 
according to Lemma |5.2j . On the other hand, z belongs to both of these sets. 
The obtained contradiction proves the theorem. □ 

6 More on Example I 

Here we describe what happens if the function a of Section 3 is slightly 
modified. Since none of the changes we consider below affects the second 
order relative spectra in any serious way, we concentrate on the standard 
(first order) projection method. 

We have performed calculations, similar to those presented in Section 
3, for the operator A corresponding to the function which is obtained from 
a by multiplying the cosines by 0.01, 0.001, or zero. The only noticeable 
difference with the results in Section 3 is that this time there are more spuri- 
ous eigenvalues in the gap of the essential spectrum and additional spurious 
eigenvalues appear, as N increases, quicker than in the original example (cf. 
Table [I] and Remark |3.1| ) . 



The situation becomes more interesting if one moves the discontinuity of 
a from to a point y ^ 0: 

i 

. v . .. . for — 7r < x < y , 

a(x) = a y {x) = { 2 2 ^ ~ y (12) 

for y < x < it . 
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The corresponding operator of multiplication by a y will be denoted A y . 

To put the results into the right perspective, we need to return to the 
original example A = A first. We know from Section 3 that Spec(A,£jv) 
show remarkable stability as N = 2n + l increases. The first surprise happens 
when one looks at Spec(v4, C N ) for N = 2n, where := span ({efc}£~i n ). 
These also seem to be converging, but the spurious eigenvalues differ consid- 



erably from those corresponding to N = 2n + 1 (see Figure [12] and Table 
and compare them with Figure |2] and Table [I]) . 



N= 1000 



Figure 12: Spec(A, £ woo ) . 



N 


100 400 700 1000 


Spec(A,£jv) n (-1,1) 


-0.9988 -0.9935 -0.9882 
-0.9910 -0.9331 -0.9033 -0.8832 
-0.7721 -0.6358 -0.5853 -0.5546 
-0.0661 0.0317 0.0636 0.0821 
0.8883 0.8539 0.8418 0.8345 



Table 2: Eigenvalues lying in the gap of the essential spectrum of A for some even 
values of N. 

One might argue that this difference between the two cases occurs due 
to the nonsymmetric truncation in the case N = 2n and that taking only 
odd values of N is more natural. We prefer, however, to interpret this fact 
as 2 -periodicity in N of the points of Spec(v4, C N ) lying in the gap (—1, 1) of 
the essential spectrum of A, i.e., of the spurious eigenvalues of A. Such an 
interpretation is supported by Table 0, which shows the spurious eigenvalues 
for some consecutive values of N. 

We show below that in a more general case of the operator A y the "period" 
depends on the location of the point y of the discontinuity of a y (see flT2])). 
This period is equal to 2 in the original example only because we take y = 
there. 

Let y = — 7r/3. Figure |13] shows the computed (spurious) eigenvalues 
of A y = a y I in the interval (—1,1) as functions of N. Remarkably, the 
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N 


301 302 303 304 305 306 307 308 


Spec(A, £ w )n(-l, 1) 


-0.9899 -0.9470 -0.9898 -0.9467 -0.9896 -0.9464 -0.9895 -0.9461 
-0.8468 -0.6622 -0.8463 -0.6616 -0.8458 -0.6610 -0.8453 -0.6603 
-0.3740 0.0143 -0.3734 0.0147 -0.3728 0.0151 -0.3721 0.0155 
0.4554 0.8604 0.4555 0.8603 0.4556 0.8601 0.4557 0.8600 



Table 3: Eigenvalues lying in the gap of the essential spectrum of A for some 
consecutive values of N. 2-periodicity in N is clearly visible. 

eigenvalues vary quite slowly with N as long as iV = const (mod 3). To 
emphasize this "3-periodicity" in N, we use different symbols for points 
with different values of N (mod 3). 
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Figure 13: Points of Spec(A y , Cn), V = — tt/3, in the interval (—1, 1) as functions 
of N. 

A similar effect occurs for the operator B y , which is obtained from the 
operator B in (||) by replacing a with a y (see Figure 0). However, curves 
corresponding to the eigenvalues in the interval (—1,1), which were confirmed 
by the second order relative spectra as the correct ones, do not depend on 
the value of N (mod 3). They also seem to converge to constants faster than 
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those depending on iV (mod 3). 
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Figure 14: Points of Spec(B y , £jv), y = — vr/3, in the interval (—1, 1) as functions 
of N. Note the presence of both "genuine" and "spurious" eigenvalues. 

We have performed similar calculations with y = —^/{p/q) = —ttq/p 
for the following values p/q = 2, 4, 5, 6, 7, 21, 21/10, 21/5, 7/3, 8/3, and 
5/2 (see Figure |I5| for the results for some of the values). In all these cases 
the spurious eigenvalues of A y showed a very clear pattern of "periodic" 
behaviour with the minimal "period" 



2, 

u(q/p) := { p, 
2p, 



if g = 0, 

if p and q are both odd, 

if q and either p or q is even. 



Such "periodic" behaviour can be observed in Figure [L5| The appearance of 
uj(q/p) here is probably not very surprising, as it is, for any simple fraction 
q/p, the smallest natural number such that ux (length of [— ir, —Trq/p}) and 
wx (length of [— nq/p, tt]) are integer multiples of 2tt. 
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Figure 15: Points of Spec(A y , Cm) in the interval (—1, 1) as functions of N, plotted 
for various values of y. 



7 More on Example II 

As we have shown in Section 4, a FEM calculation of the spectrum of 
operator H using the projection space L N x h N leads to heavy spectral 
pollution in the interval (1,2). At the same time, the second order spectra 
relative to the projection space M^r x L^v approximate the spectrum of H 
very well. 

A natural question is whether such a discrepancy may be mostly due 
to the different projection spaces being used in the calculations of the first 
order and second order spectra. To check this, we compute the spectra 
Spec(if, M^v x Ln) using the same projection spaces as we used in Section 
4 for the second order relative spectra. The results (restricted to eigenvalues 
lying in the interval [0,6]) are shown in Figure [TJ]. This time, there are no 
spurious eigenvalues at all! In the particular case N = 100, we have 100 
eigenvalues in [0,1], and 200 eigenvalues in [2, oo) (recall that din^M^ x 
Ljv) = 3N). The actual eigenvalues of H are well approximated in this case 
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as well, see Remark 17.11 below 
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Figure 16: FEM calculated Spec(iT, x Ljv) plotted against N. There are no 
spurious eigenvalues! 

This of course indicates that the projection spaces Mjv x are much 
better suited for the calculation of the spectrum of H than the projection 
spaces Ltv x L/v- Such a dependence upon the choice of the spaces of finite 
elements is quite interesting, especially because it does not occur for the 
operator —d 2 /dx 2 and because of the explicit relations which exist between 
the eigenvalues and eigenfunctions of H and those of —d 2 /dx 2 , see (§) and 

The spaces Mjy x Ljv are well suited for the calculation of Spec(-ff) be- 
cause they are subspaces of the mixed order Sobolev space 1F 2 2 ([— n, ir))) x 
PF 2 1 ([— 7r, 7r)]) and therefore reflect much better the mixed orders of differ- 
entiation in the operator H, see also [[RSHSPVf and references therein. We 



attempt to explain this phenomenon in a more general setting below, see 
Remark 



One may be tempted to think that the success of the second order method 
used in Section 4 is also only due to the correct choice of the finite elements. 
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In order to check this assumption let us compare the results of calculation of 
the first order and second order relative spectra using the projection spaces 
Mjy x Mjy. The results for Spec(if, M^r x M^v) are shown in Figure |T7 - 
again, there is spectral pollution in the interval (1,2). 



1 i ; i i i i i i i M i i i i i i i i i i i i i i i 



10 15 20 25 30 35 40 45 50 55 60 



Figure 17: Spec(H, M^v X Mjv) plotted against N. Again, all the eigenvalues in 
the interval (1,2) are spurious! 

On the other hand, the calculations of the second order relative spectra 



Spec 2 (-£f, Mat x M^) produce pictures very similar to those in Figures |10 
and [D]. The only noticeable difference is the presence of a number of points 
with relatively large moduli of imaginary parts (> 10 for iV = 50) which seem 



to lie on a "nice" curve passing through the point A = 1 (see Figure y| note 
the different scale along the imaginary axis compared with either Figure |10 
or 



11) 



Thus, for the spaces x Mjv, the projection method using second order 
relative spectra works quite well, whereas the standard method again fails 
dramatically. 
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Figure 18: Spec 2 (#, M 50 x M 50 ). 

Remark 7.1. It should be noted that both the first order and the second or- 
der relative spectra (outside (1,2)) approximate Spec(iJ) very well, with ac- 
curacy depending on the finite elements chosen. The points of Spec(if, L^r x 
Ljv), Spec(if, Mjv x L N ), and Spec(if, M^r x M^r) converge to the points 
A^ E Spec(iJ) with the rate of about 0(N- 2 ), 0(iV- 4 ), and 0(N' 6 ), respec- 
tively For the second and third choice of finite elements, the same rates of 
convergence are achieved by the real parts of the points of Spec 2 (if, xLjv) 
and Spec 2 (if, MjvxMj/), whereas the imaginary parts converge to zero with 
the rate of about 0(N~ 2 ) and 0(N~ 3 ), respectively. These convergence rates 
were obtained by analyzing the log-log dependence of errors upon TV for sev- 
eral low eigenvalues taken from both A + and A - series. 

Remark 7.2. Here we try to explain why the standard projection method 
of approximating Spec (if) does not lead to spectral pollution in the case of 
the projection spaces x Ljy. The main reason is that it seems, in this 
particular case, to "mimic" the way one finds Spec (if) analytically. Indeed, 
let C\ be a space of finite elements for v , C 2 be a space of finite elements 
for u, and let £J 2 '■= {if)' '■ if) E C 2 }. Let ipj, j = be a basis of Ci 
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and ipki k — 1, . . . , N be a basis of £2- Since the projection method does not 
depend on the choice of bases, provided all the calculations are precise, we 
can assume that the basis {(fij} is orthonormal. Then the finite dimensional 
system for (|6|) takes the form 



m = l,...,N, 



^ 1 (v./- <7,J''/, + -<'., 

./ 1 Vf. 

This leads to the following system for the coefficients u^. 



N N M 



k=l k=l j=l 

N 

A y^X^m, ^k)Uk, m = l,...,N. (13) 
k=l 

Let Pi be the orthogonal projection onto C\. Then 

M M 

Pl^'k = and (PAi Pl^'k) = ^2Wmi VMPs' tic)- 

3=1 3=1 

Hence ( |T3"D is equivalent to 

N N N 

X^(V'm) ^'k) U k + ^(Pl^L PM^k = A ^{ipm, IpkW, (14) 



A-2 

fc=l fc=l fc=l 



m = 1 JV. 



On the other hand, it is clear that (0) implies 



-u"-^u" = \u, (15) 



which is equivalent to (|7|). So, the projection method applied to ([15]) should 
not lead to spectral pollution. This gives the following finite dimensional 

system 

N I N N 

X^m' ^'k) U k + V>[K = A ^(VVn, ^fcK, (16) 



A-2 

k=l k=l k=l 



m = 1 iV. 



If C\ ~D C' 2 -, then (|T4j) and ( |16|) coincide and the projection method applied 
to (|6|) does not pollute. If C\ 2 A, Du ^ £1 an d £' 2 are "compatible", then 
it may well happen that (|14D mimics the properties of (|16D at least partially, 
and the projection method applied to (ffi) does not lead to spectral pollution. 
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8 Concluding remarks 



The main message of the paper is that projection methods using second 
order relative spectra may provide an efficient tool for dealing with spectral 
pollution. They are not much harder to implement than the standard pro- 
jection methods and are likely to give useful a posteriori information (see 
Theorem |2.5| ) about the location of (at least a part of) the spectrum of a 
self-adjoint operator. We would feel that the paper has achieved its goal if 
it persuaded other people to try second order relative spectra in their prob- 
lems of interest. We hope that the method we suggest here will work well in 
problems more serious than the model examples considered above. 

We would like to stress again that we do not claim that second order 
relative spectra are the only way of addressing spectral pollution. In our 
Example I, for instance, one could argue that spectral pollution could have 
been detected just by looking at pictures like Figure 3, without any reference 
to the second order relative spectra. The argument would be that a genuine 
eigenfunction "cannot" look like the one in Figure 3, and hence the corre- 
sponding eigenvalue is spurious. The "periodicity" of spurious eigenvalues 
(see Section 6) might also have been used to detect possible spectral pollu- 
tion. Indeed, those approximate eigenvalues of B y , which were confirmed by 
the second order relative spectra as the correct ones, do not exhibit "period- 
icity" in N. They are present for any sufficiently large value of N. We have 
no intention to disagree with these arguments. Note however that both of 
these approaches are based on the results of calculations using the standard 
(first order) projection method, which does not see the difference between A 
and the Toeplitz operator T a . Consequently, they cannot distinguish between 
spectral pollution for the former and non-pollution for the latter. Second or- 
der relative spectra, on the other hand, are free from this problem. Our point 
is that it is very easy to fail to notice spectral pollution unless one is very 
careful, and that second order relative spectra seem to be very efficient "lie 
detectors" . 

Another perfectly legitimate argument is that one should try to avoid 
spectral pollution by making a sensible choice of the subspaces rather 
than allow it to happen and then try to detect it. The standard projection 
method using x L/v, for instance, works very well in our Example II 
(see Section 7). One should, however, bear in mind that choosing "right" 
subspaces when applying the standard projection method in more difficult 
situations may prove non-trivial and require a good understanding of the 
problem at hand. Second order relative spectra, on the other hand, seem to 
provide a more robust projection method which works well even in the case 
of "wrong" subspaces, when the standard one fails. The case of M^r x M. N 
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considered in Section 7 is an example of such situation. 

Finally, we refer the reader to the companion paper by E. B. Davies and 
M. Plum |D~F ], where a different method of treating spectral pollution is 
developed. 
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